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Abstract 

We revisit the classic problem of the structure of the isotropic-nematic interface within Ginzburg- 
Landau-de Gennes theory, refining previous analytic treatments of biaxiality at the interface. We 
compare our analysis with numerical results obtained through a highly accurate spectral collocation 
scheme for the solution of the Landau-Ginzburg-de Gennes equations. In comparison to earlier 
work, we obtain improved agreement with numerics for both the uniaxial and biaxial profiles, 
accurate asymptotic results for the decay of biaxial order on both nematic and isotropic sides of 
the interface and accurate fits to data from density functional approaches to this problem. 

PACS numbers: 42.70.Df,67.30.hp,61.30.Dk,61.30.Hn 
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Liquid crystalline states of matter provide a useful testing ground for statistical me- 
chanical theories of interface structure, since a variety of ordered phases can be accessed in 
experiments and computer simulations. The structure of the isotropic-nematic (I-N) inter- 
face presents a simple example of how interfacial order can differ radically from order in the 
coexisting bulk phases, since biaxial order is generically expected at the interface even if the 
stable ordered phase is purely uniaxial. The study of the isotropic-nematic interface was 
initiated in an insightful paper by de Gennes, who introduced a simple uniaxial ansatz for 
the tensor order parameter Qap which describes nematic order Ij. The de Gennes ansatz is 
exact in the absence of elastic anisotropy. However, the description of the interface in the 
presence of such anisotropy poses a formidable analytic and numerical problem, since the 
partial differential equations for the five independent components of Qai3 contain non-linear 
couplings, while Qajj is itself constrained by symmetry and the requirement that its trace 
vanish. 

Popa-Nita, Sluckin and Wheeler (PSW) ^ studied the I-N interface incorporating elastic 
anisotropy in the limit of planar anchoring, adapting a parametrization introduced by Sen 
and Sullivan[3]. In this parametrization, the principal axes of Qaf3 remain fixed in space, and 
the problem reduces to the solution of two coupled non-linear partial differential equations 
in the dimension perpendicular to the interface. These equations represent the variation 
of the amplitude of uniaxial and biaxial ordering across the interface. PSW showed that 
the solutions of these equations exhibited biaxiality in a region about the interface The 
uniaxial order parameter (S) was adequately represented by a tanh profile, as in the original 
calculation of de Gennes, while the biaxial order parameter (T) exhibited more complex 
behaviour, peaking towards the isotropic side and with a trough on the nematic side. The 
biaxial profile was also shown to have a long tail towards the isotropic side, a feature hard 
to anticipate on physical grounds. 

This paper extends these calculations in several new ways. First, we show that terms 
dropped by PSW in their simplification of the Ginzburg-Landau-de Gennes (GLdG) equa- 
tions are, in fact, comparable in magnitude to the terms they retain, especially for small 
values of k = L2/L1, the ratio of the coefficients of the two lowest-order gradient terms in 
the GLdG expansion. Thus, a more accurate treatment of the interface requires that these 
terms be retained. The resulting equations have closed form solutions in terms of hypergeo- 
metric functions. We show that such solutions provide a better description of the numerical 



2 



data than the original calculation of PSW. We benchmark our analytic results through an 
accurate numerical procedure, based on a Chebyshev polynomial expansion, for the study 
of these equations. 

We begin with the GLdG expansion of the free energy for a general Q^js 

+ ^Li{daQl3'r){daQl3'y) + ^L2{daQa'r){d[^Q f^^)]. (l) 

Here A,B and C are expansion parameters, while Li,L2 are elastic constants. We choose 
B = —0.5, C = 2.67 and A = _B^/27C, thus enforcing phase coexistence between an isotropic 
and uniaxial nematic phase IJ] . The interface is taken to be fiat and infinitely extended in the 
x — y plane. The spatial variation of the order parameter only occurs along the z directionjsl. 
We scale Qa/3 Qap/Sc where S'c = — ^ 9^%-^' ^"^^ measure lengths in units of 

Ic = v^54C(Li + 2L2/3)/52; we choose Li = IQ-^ in our numerics and obtain L2 from our 
choice of k. In the case of planar anchoring, the ordering at infinity is purely uniaxial and 
taken to be along the x axis. In this shown by Sen and Sullivan, uniaxial and biaxial 

order vary only with z and the principal axes of the Q tensor remain fixed in space. The 
form of Q is then 

^ 

Q= |(-5 + T) 

^0 -1{S + T)J 

Inserting this form of Q into the free energy and performing the minimization yields 

(6 + '').92s + TTrArT^'^ = 45 - 125' + + 4T^ + (3) 



(2) 



(3 + 2k) ' {3 + 2k) 



,2c , (2 + «^) ^2^ _ 4^ , , 8T3 , 8S'T 



-diS + ; \ dtT = -T + 8ST + + . (4) 

(3 + 2k) ^ (3 + 2k) ^ 3 9 3 ^ ^ 

Popa-Nita, Sluckin and Wheeler now make several approximations to Eqs. [3] and |4] to 
solve them. First, in Eq. [3l all terms in T are dropped, since S is typically much larger than 
T. The resulting equation for 5* is solved by the tanh function. In Eq. HI PSW drop the 
(3+2 k) '^z'^ term while retaining (2,+2k) ^'z^- ^ ^^^^ self-consistency of this approximation is 
the comparison of the magnitude of these terms within the theory. Fig. [T] shows the terms 
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FIG. 1: (Color online) A comparison of the terms -^i^^d'^T (dark line) and j^^j^d'^S (light 
line) obtained within the PSW solution to the GLdG equations, for a k value of 18.0. The PSW 
approximation consists of ignoring the j§^^d'^T term in comparison to the ^gq^^^^S* term. Both 
terms, however, are of comparable magnitude. 

{3+2k) (dark line) and j^^2k)'^z'^ (light line) computed through the PSW solution. As 
can be seen from the figure these terms only differ by a factor of order unity. Deep into 
the isotropic side, the term ignored by PSW exceeds the value of the term retained. Thus, 
while the PSW approach leads to a straightforward algebraic relation between T and S, a 
more accurate method would be to retain the partial derivative term as well, requiring that 
we solve a partial differential equation as opposed to an algebraic one. 

Our approach to this problem uses the same approximations as PSW for Eq. [3l We thus 
take 



Q _ 



1 + tanh( 



(5) 



where ^ 



Inserting this in equation (jlj), scaling z by v^^, redefining the resulting 



1+2k/3 

quantity as z again, and dropping the nonlinear term, we obtain 



dlT = 2/?[tanh^(2) + 8tanh(2) + 9]T + 



2 + K 



tanh{z) [1 + tanh{z)] [1 - tanh(z)] . (6) 



with j3 



6+K 

3(2+k) ■ 




-1 -0.5 0.5 1 -0.2 -0.1 0.1 




FIG. 2: (Color online) Biaxial and uniaxial profiles for k = 0(a), 0.4(6), 4(c) and 18.0(d), comparing 
results from our numerical computations (x), with our analytic formula (dashed line) and the 
formula of PSW (solid line). The main figure shows the biaxial profile whereas the inset shows 
the uniaxial profile. In (a), for k = 0, the solution has T = 0, with the S profile exactly given by 
the tanh form. In (b), for k = 0.4, the computed biaxial profile (T) (main panel) is fit remarkably 
well by our analytic form, whereas the PSW approximation tends to overestimate the peak value. 
The uniaxial (S) profile is shown in the inset of (a); here the results obtained by us and by PSW 
are identical and the fit to a tanh profile is accurate over the entire region. In (c) (main panel), 
for K = 4.0, the numerical data are fit well by the analytic forms, particularly away from the main 
peak, yielding essentially exact agreement deep into the isotropic and nematic sides. The PSW 
approximation is still an overestimate to the peak value, and also differs sharply in relation to the 
numerical data deep into the isotropic side. The inset shows the uniaxial (S) profile for this case. 
In (d) (main panel), for k = 18.0, the PSW form appears to fit the peak better for larger k, but 
again fails to capture the decay towards the isotropic side. 
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FIG. 3: (Color online) A comparison of the results of our analytic calculation to profiles of T 
obtained from a density functional calculation for the isotropic-nematic interface. Profiles obtained 
for two values of k, k = 5.8 (for z < 0) and k = 0.69 (for z > 0) are shown. The larger k value 
essentially fits the T profile exactly on the isotropic side, whereas the smaller k value provides an 
accurate fit on the nematic side. The inset shows the S profile obtained from the density functional 
calculation, together with an optimum fit varying the value of Ic 
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The solution of the equation consists of a homogeneous part Th and a particular part 
Tp where n = C^t^^il - ty^^Fila^, h, Ci, t] + C2t-^^{1 - t)^2Fi[a2, h, C2, 1 - t] and 

i-2v^+^^; C2 = l + 2v/;5andt= ( 
function and Ci and C2 are fixed by boundary conditions. 
The particular solution takes the form 



ci = 1 + 6v^;Q2 = 9~^v/^~ '9 1 ^2 
The function 2F1 is a hypergeometric 



with 



T,{z) = [-yi{z)l2{z) + y2{z)h{z)] /W{z), 



1 -tanh(2)\3v^/l + tanh(z)\-v^ . , 1 - tanh(z) 
( 7, ) 2i^i[ai,0i,ci, 



(7) 



2 / V 2 
1 - tanh(2) -3v^ / 1 + tanh(;z) >^ 

\ 2 



Fi[a2,&2,C2: 



1 + tanh(2;) 



(8) 
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and the Wronskian W{z) = W = yi{dy2/dz) — y2{dyi/dz), where 
2/t \ - (oi) 

m 



2 + K ^ (ci)„m! V 2 + m + 2^ 

m=0 



m+^VP p [1 2 + m + 2Jp, 2 + m + 3 v^, tl) 



(9) 



h(z) 



2k \- (Q.2)n(^2)n ^l+n+v^^-| ^ ^l-3^/:a 
n=0 



2 + /^^^^ (C2)nr2! ^ ^ V 2 + n-2v^ 

^ 2F1 [1, 2 + n - 2 v^, 2 + n + v<5, 



(2 + n-2v^)(l + n + v^)- 



(10) 



The Pochhammer symbol (a)„ which enters above is defined via (a)„ = a(a+l)(a+2) . . . (a+ 
n— 1). Here ti = [1 + tanh(2;)]/2 and the result for Ii{z) and hiz) is obtained by expanding 
the hypergeometric functions in Eqns. [H]in a power series and integrating term-by-term [sl- 
Note that the solutions of the homogeneous part diverge asymptotically. Thus, for the 
boundary condition T = at z = ±00 the only physical solution is the particular one. Eq. [7] 
is thus the key analytical result of this paper, describing the variation of biaxiality across 
the interface. In our numerical evaluations, we sum the series for Ii{z) and l2{z), retaining 
as many terms as are required to ensure convergence. The series in I2 converges very fast 
(only 3 terms need be retained for good results) whereas the series in Ii converges more 
slowly and around 9 terms must be retained for convergence. To convert these into physical 
units, we must undo the sequence of length transformations, replacing z zj[^2^l^. 
An asymptotic analysis of these equations is possible: for z —00, 5* and T are small. 

2z 

The tanh profile for S can be approximated as |(1 -|-tanh(^^)) ^ e^^? while Eq. (jl]) takes 
the form 2edlT = 4/?T- (^)e^^ with P = Thus d^T = I^T- (£|fgye^^ 

with asymptotic solution 

r ^ eV 3(2+.) _^ _oo, (11) 

a result in perfect accord with the computed forms of T deep into the isotropic phase. Note 
that T > S and (S + T)/2 > S as one moves deeper into the isotropic side. This implies 
that the principal order parameter is negative as pointed out in Ref. where this result 
was obtained numerically. 



7 



As z ^ oo, an alternative asymptotic expansion can be derived by taking S* = 1 — |e 
with b = 2,/4r^. We then obtain 



_9 /3+2k_ 

T ~ e 6+-^,2 oo, (12) 

in agreement with our numerical results. Popa-Nita, Sluckin and Wheeler provide an analysis 
of the asymptotics in the specific limit that /t ^ oo. However, our results cannot be directly 
translated to this limit, since we assume a tanh profile of S; this approximation becomes 
increasingly inaccurate for larger k (see below). 

Our numerical results are obtained using a spectral collocation method [6] , applied to our 
knowledge for the first time to the GLdG equations. In the spectral collocation, the solution 
is expanded in an orthogonal basis of Chebyshev polynomials in a bounded interval. Dif- 
ferentiation operators constructed from this Chebyshev interpolant are spectrally accurate, 
in the sense that the error vanishes exponentially in the number of retained polynomials. 
The interpolant is constructed so as to satisfy Dirichlet boundary conditions. Though the 
physical problem is for an unbounded interval, our numerical approximation of a bounded 
interval gives excellent results since all variation in the order parameters is restricted to the 
region proximate to the interface. 

Specifically, we solve the equations of equilibrium 



{A + CTrQ^)Q^p{^,t) + BQlp{^,t) = LiV2g«/j(x, t) + L2 V„(V^g/37(x, t)) (13) 
by transforming to a basis {oj} which enforces symmetry and tracelessness, as Qa/3 = 



ELi«i^a/3. where, = y'|zz,T2 = ^l(xx-yy),T3 = v^xy,T4 = v^xz, 
= y z . Overbars indicate traceless symmetric parts. We thus obtain five simulta- 
neous partial differential equations for the a^, which are steady-states of the time- dependent 
equations we have obtained earlier 3]- Note specifically that we make no symmetry-based 
ansatz for the components of Qap \a\. 

The spectral collocation reduces these differential equations to non-linear algebraic equa- 
tions. We solve them using a relaxation method from a well-chosen initial condition, relaxing 
till the differential change in successive iterations is less than 10~^. Spectral convergence to 
machine accuracy is obtained by retaining 128 Chebyshev modes, as we have checked by an 
explicit calculation. To compare with analytical and density functional results, the solution 
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at the Chebyshev nodes is interpolated using barycentric interpolation without compromis- 
ing spectral accuracy. The DMSUITE library is used for the numerical implementation 

Our results are summarized in Fig. [2] and Fig. [31 The main panel of Fig. [2] (a), obtained 
by solving Eq. [3] and H] for a value of k = 0.0, shows the biaxiality profile obtained using 
our numerical spectral scheme (crosses), as compared to the analytic result of T = 0. The 
uniaxiality profile shown in the inset is exactly the tanh profile obtained by de Gennes. This 
limit provides a simple test of our numerical methods, since the solution to Eqs. [3] and H] 
in this limit is exact. Fig. [2] (b) shows the biaxiality profile obtained using our spectral 
scheme (crosses), as compared to the analytic results derived here (dashed line) and results 
obtained by PSW (solid line) for a value of /t = 0.4. As can be seen, the numerical data 
are fit remarkably well by the analytic forms, whereas the PSW approximation tends to 
overestimate the peak value. The inset to Fig. [2] (b) shows the uniaxial (S) profile, obtained 
numerically as well as in our analytic calculation; here the results obtained by us and by 
PSW are identical. The fit to a tanh profile is accurate over the entire region. 

The main panel of Fig. [2] (c) shows the biaxiality profile obtained using our spectral 
scheme (crosses), as compared to the analytic results derived here (dashed line) and results 
obtained by PSW (solid line) for a value of k = 4. Again the numerical data are fit well by the 
analytic forms, particularly away from the main peak, yielding essentially exact agreement 
deep into the isotropic and nematic sides. The PSW approximation is still an overestimate 
to the peak value, and also differs sharply in relation to the numerical data deep into the 
isotropic side. The inset to Fig. [2] (c) shows the uniaxial (S) profile for this case. Fig. [2] (d) 
shows the biaxiality profile obtained using our spectral scheme (crosses), as compared to the 
analytic results derived here (dashed line) and results obtained by PSW (solid line) for a 
value of K = 18. For these - and larger - values of k, our analytic fits differ noticeably from 
the numerical data. The PSW form appears to fit better for larger k, although we believe 
that this is fortuitous. It appears that the principal error arises from our approximation of 
the S profile as a tanh form. For large k, this approximation is less accurate. 

Fig. [3] compares the results of our analytic calculation to profiles of T obtained from a 
density functional calculation for the isotropic-nematic interface 10(] a method which pro- 



vides an alternative, more molecular approach to this problem 



We have taken numerical 



data for uniaxial and biaxial profiles obtained in Ref. lOj, varying the free parameters Sc, 
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Ic and K in our solutions till an optimal fit is obtained. The values of Sc and Ic can be 
obtained from fits to 5*; thus only k need be varied to represent the T profile. Fig. [3] 
shows profiles obtained for two values of k: k, = 5.8 (for z < 0) and k = 0.69 (for z > 0). 
The larger k value fits the profile very closely on the isotropic side, whereas the smaller k 
value provides an accurate fit on the nematic side [l^. It does not seem possible to fit the 
complete profile using a single value of k. This could have been anticipated on physical 
grounds since the density functional theory yields a density difference between coexisting 
isotropic and nematic phases. The elastic coefficients Li and L2 which enter our calculation 
do in principle contain a density dependence which we ignore here. 

In conclusion, we have presented results for the uniaxial and biaxial profiles, in the case of 
planar anchoring, for the classic problem of the structure of the isotropic-nematic interface 
within Ginzburg-Landau-de Gennes theory. Our work refines previous analytic treatments 
of biaxiality at the interface. We have implemented a highly accurate spectral collocation 
scheme for the solution of the Landau- Ginzburg-de Gennes equations and used this numerical 
scheme in our tests of the analytic results. 

In comparison to earlier work, we obtain improved agreement with numerics for both the 
uniaxial and biaxial profiles, with our results being increasingly accurate as the anisotropy 
is reduced. We also provide accurate asymptotic results for the decay of the S and T 
order parameters deep into the nematic and isotropic phases. Our calculated profiles show 
a pleasing consistency with profiles obtained from density functional approaches. Further 
extensions of these numerical and analytic methods to the case of an intermediate anchoring 
condition far from the interface are currently under way. 
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the DST (India) and the Indo- French Centre for the Promotion of Advanced Research. 
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